### Replication code for Taylor C. Boas, F. Daniel Hidalgo, and Guillermo Toral. "Competence versus Priorities: Negative Electoral Responses to Education Quality in Brazil"
### This function helps visualize regression discontinuity designs
### R version, platform, and package versions reported at the end of the file
### June 19, 2020    

# rd_plot: FUNCTION TO MAKE RD PLOTS --------------------------------------

rd_plot <- function(x,y,cutpoint=0,maintitle="",xtitle="",ytitle="",xrange=c(-1,1),
                     yrange=c(0,1),bins=100,legendabove="",legendbelow="", titlesize=1){
  lobelow <- loess(y~x,subset = x < cutpoint)
  loabove <- loess(y~x,subset = x > cutpoint)
  xl <- seq(min(x,na.rm=T),max(x,na.rm=T), (max(x,na.rm=T) - min(x,na.rm=T))/1000)
  xbelow <- seq(xrange[1],cutpoint,(xrange[2] - xrange[1])/900)
  xabove <- seq(cutpoint,xrange[2],(xrange[2] - xrange[1])/900)
  pred.below <- predict(lobelow,xbelow,se=T)
  pred.above <- predict(loabove,xabove,se=T)
  minbelow <- pred.below$fit - pred.below$s*1.96
  maxbelow <- pred.below$fit + pred.below$s*1.96
  minabove <- pred.above$fit - pred.above$s*1.96
  maxabove <- pred.above$fit + pred.above$s*1.96
  xsubset <- subset(x,x<xrange[2] & x > xrange[1])
  sub <- as.data.frame(cbind(x,y))
  sublow <- sub[which(sub$x > xrange[1] & sub$x < cutpoint),]
  subhigh <- sub[which(sub$x < xrange[2] & sub$x > cutpoint),]
  require(Hmisc)
  minx <- rep(NA,bins)
  maxx <- rep(NA,bins)
  meany <- rep(NA,bins)
  means <- as.data.frame(cbind(minx,maxx,meany))
  for(i in 1:(bins/2)){
    means$minx[i] <- min(split(sublow, cut2(sublow$x, g=bins/2))[[i]][[1]],na.rm=T)
    means$maxx[i] <- max(split(sublow, cut2(sublow$x, g=bins/2))[[i]][[1]],na.rm=T)
    means$meany[i] <- mean(split(sublow,cut2(sublow$x, g=bins/2))[[i]][[2]],na.rm=T)
  }
  for(i in 1:(bins/2)){
    means$minx[i+bins/2] <- min(split(subhigh, cut2(subhigh$x, g=bins/2))[[i]][[1]],na.rm=T)
    means$maxx[i+bins/2] <- max(split(subhigh, cut2(subhigh$x, g=bins/2))[[i]][[1]],na.rm=T)
    means$meany[i+bins/2] <- mean(split(subhigh,cut2(subhigh$x, g=bins/2))[[i]][[2]],na.rm=T)
  }
  means$center <- (means$minx+means$maxx)/2
  plot(x,y,pch=20,col="lightgray",main=maintitle,
       xlab=xtitle,
       ylab=ytitle,type="p",cex=.8,
       xlim=xrange,ylim=yrange, cex.lab=1.2, cex.main=titlesize)
  rug(x,col="black",ticksize=0.03,lwd=0.3)
  polygon(c(xabove,rev(xabove)),c(maxabove,rev(minabove)),
          col=adjustcolor("darkgreen",alpha.f=0.3),border="white")
  polygon(c(rev(xbelow),xbelow),c(rev(maxbelow),minbelow),
          col=adjustcolor("red",alpha.f=0.3),border="white")
  abline(v=cutpoint,col="darkblue",lty=2,lwd=2)
  lines(means$center[1:bins/2],means$meany[1:bins/2],col="red",type="p",pch=19,cex=.8)
  lines(means$center[bins/2+1:bins],means$meany[bins/2+1:bins],col="darkgreen",type="p",pch=19,cex=.8)
  lines(xbelow, pred.below$fit, col='black', lwd=2)
  lines(xabove, pred.above$fit, col='black', lwd=2)
  legend(x=cutpoint,y=yrange[2],legend = legendabove, bty ="n", pch=NA, text.col="darkgreen",text.width=3, cex=1.2) 
  legend(x=xrange[1],y=yrange[2],legend = legendbelow, bty ="n", pch=NA, text.col="red",text.width=3, cex=1.2) 
}

# NOTES -- R version, platform, and loaded packages -------------------------
# sessionInfo(package = NULL)
# R version 3.6.3 (2020-02-29)
# Platform: x86_64-apple-darwin15.6.0 (64-bit)
# Running under: macOS Catalina 10.15.3
# 
# Matrix products: default
# BLAS:   /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib
# LAPACK: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRlapack.dylib
# 
# Random number generation:
#   RNG:     Mersenne-Twister 
# Normal:  Inversion 
# Sample:  Rounding 
# 
# locale:
#   [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
# 
# attached base packages:
#   [1] stats     graphics  grDevices utils     datasets  methods   base     
# 
# other attached packages:
#   [1] Hmisc_4.4-0       lattice_0.20-40   timelineS_0.1.1   lfe_2.8-5         Matrix_1.2-18    
# [6] xtable_1.8-4      texreg_1.37.1     rdrobust_0.99.7   rdd_0.57          Formula_1.2-3    
# [11] AER_1.2-9         survival_3.1-11   car_3.0-7         carData_3.0-3     lmtest_0.9-37    
# [16] zoo_1.8-7         sandwich_2.5-1    shiny_1.4.0.2     readxl_1.3.1      codebook_0.9.2   
# [21] electionsBR_0.3.1 forcats_0.5.0     stringr_1.4.0     dplyr_1.0.0       purrr_0.3.3      
# [26] readr_1.3.1       tidyr_1.0.2       tibble_3.0.0      ggplot2_3.3.0     tidyverse_1.3.0  
# 
# loaded via a namespace (and not attached):
#   [1] colorspace_1.4-1    ellipsis_0.3.0      rio_0.5.16          htmlTable_1.13.3   
# [5] base64enc_0.1-3     fs_1.4.1            rstudioapi_0.11     listenv_0.8.0      
# [9] farver_2.0.3        DT_0.13             fansi_0.4.1         lubridate_1.7.4    
# [13] xml2_1.3.0          codetools_0.2-16    splines_3.6.3       knitr_1.28         
# [17] jsonlite_1.6.1      broom_0.5.5         cluster_2.1.0       dbplyr_1.4.2       
# [21] png_0.1-7           compiler_3.6.3      httr_1.4.1          backports_1.1.5    
# [25] assertthat_0.2.1    fastmap_1.0.1       cli_2.0.2           later_1.0.0        
# [29] acepack_1.4.1       htmltools_0.4.0     tools_3.6.3         gtable_0.3.0       
# [33] glue_1.3.2          Rcpp_1.0.4          cellranger_1.1.0    vctrs_0.3.1        
# [37] nlme_3.1-145        crosstalk_1.1.0.1   xfun_0.12           globals_0.12.5     
# [41] openxlsx_4.1.4      rvest_0.3.5         mime_0.9            miniUI_0.1.1.1     
# [45] lifecycle_0.2.0     future_1.17.0       scales_1.1.0        hms_0.5.3          
# [49] promises_1.1.0      parallel_3.6.3      RColorBrewer_1.1-2  yaml_2.2.1         
# [53] curl_4.3            gridExtra_2.3       rpart_4.1-15        labelled_2.5.0     
# [57] latticeExtra_0.6-29 stringi_1.4.6       highr_0.8           checkmate_2.0.0    
# [61] rmdpartials_0.5.8   zip_2.0.4           repr_1.1.0          rlang_0.4.6        
# [65] pkgconfig_2.0.3     evaluate_0.14       htmlwidgets_1.5.1   labeling_0.3       
# [69] tidyselect_1.1.0    magrittr_1.5        R6_2.4.1            generics_0.0.2     
# [73] DBI_1.1.0           pillar_1.4.3        haven_2.3.1         foreign_0.8-76     
# [77] withr_2.1.2         nnet_7.3-13         abind_1.4-5         modelr_0.1.6       
# [81] crayon_1.3.4        utf8_1.1.4          rmarkdown_2.1       jpeg_0.1-8.1       
# [85] grid_3.6.3          data.table_1.12.8   reprex_0.3.0        digest_0.6.25      
# [89] httpuv_1.5.2        munsell_0.5.0       skimr_2.1.1      